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Abstract 

Background: In a rapidly changing world, it is of fundamental importance to understand processes constraining or 
facilitating adaptation through microevolution. As different traits of an organism covary, genetic correlations are expected 
to affect evolutionary trajectories. However, only limited empirical data are available. 

Methodology/Principal Findings:YJe investigate the extent to which multivariate constraints affect the rate of adaptation, 
focusing on four morphological traits often shown to harbour large amounts of genetic variance and considered to be 
subject to limited evolutionary constraints. Our data set includes unique long-term data for seven bird species and a total of 
10 populations. We estimate population-specific matrices of genetic correlations and multivariate selection coefficients to 
predict evolutionary responses to selection. Using Bayesian methods that facilitate the propagation of errors in estimates, 
we compare (1) the rate of adaptation based on predicted response to selection when including genetic correlations with 
predictions from models where these genetic correlations were set to zero and (2) the multivariate evolvability in the 
direction of current selection to the average evolvability in random directions of the phenotypic space. We show that 
genetic correlations on average decrease the predicted rate of adaptation by 28%. Multivariate evolvability in the direction 
of current selection was systematically lower than average evolvability in random directions of space. These significant 
reductions in the rate of adaptation and reduced evolvability were due to a general nonalignment of selection and genetic 
variance, notably orthogonality of directional selection with the size axis along which most (60%) of the genetic variance is 
found. 

Conclusions: These results suggest that genetic correlations can impose significant constraints on the evolution of avian 
morphology in wild populations. This could have important impacts on evolutionary dynamics and hence population 
persistence in the face of rapid environmental change. 
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Introduction 

With the realisation that evolution can occur rapidly, there has 
been growing interest in measuring short term microevolutionary 
responses to natural selection and attempting to predict such 
responses in wild populations [1]. Understanding responses to 
selection is an exciting challenge, both at the fundamental level 
when attempting to understand and predict evolutionary mech- 
anisms and at the applied level, as in the case of management of 
responses to anthropogenic changes such as global warming [2]. 



However, the relevance of our predictions for evolutionary 
trajectories in natural settings will depend on how accurately we 
can assess selective pressures and evolutionary potential. 

Evolutionary potential is often estimated as heritability (h 2 ). 
However, most studies have reported a discrepancy between 
predicted and observed evolutionary responses to selection when 
using the breeder's equation [3,4], where heritability is multiplied 
by selection to obtain the expected evolutionary response [5] . One 
of the possible explanations for this discrepancy is that the 
estimates of evolutionary potential are inaccurate [6]. A major 
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Figure 1. Measures of constraints on response to selection for two traits, z n and z 2 . G is represented by the ellipse, g max is the first 
eigenvector of G, p is the vector of directional selection, and Az g is the response to selection calculated from the multivariate breeder's equation in 
the presence of genetic correlations. e p , the multivariate evolvability, is the projection of the response to selection on (3. e grnax is the angle between 
g max and p. Redrawn from [12]. 
doi:10.1371/journal.pone.0090444.g001 



limitation of equating the evolutionary potential of a character 
with its heritability comes from the fact that phenotypes result 
from the interaction of several characters that are functionally, 
developmentally and genetically linked. Approaching phenotypes 
as a set of independent traits may thus give a very misleading 
picture of expected phenotypic responses to selection [7,8]. Hence, 
estimating evolutionary potential requires understanding and 
assessing how genetic architecture, notably genetic correlations 
between traits, influences responses to selection, either by 
constraining or by facilitating such responses [9, 1 0] . 

Technically, such an assessment implies the estimation of 
genetic correlations as well as selection on correlated characters. 
The G matrix, the matrix of additive genetic variances and 
covariances, summarizes the genetic architecture for a set of traits. 
A geometrical representation of the G matrices can help to 
visualize how the information contained in the G matrices can be 
interpreted in terms of evolutionary potential for a given set of 
traits in a population. A spherical G matrix (with equal amount of 
additive genetic variance in all directions) provides an opportunity 
for the same amount of evolutionary response in all directions of 
phenotypic space. An elliptical G matrix, on the other hand, is 
characterized by a main axis of additive genetic variance (g max , 
Fig. 1, [11]). This axis represents the direction of highest 
evolvability in the phenotypic space, hence the direction in which 
an evolutionary response is facilitated [12]. If selection and 
are not aligned, the response to selection will be slower and 



evolvability reduced. Multivariate evolutionary constraints arise 
when there is little or no genetic variation in the direction of 
selection i.e. if some traits are genetically negatively correlated but 
submitted to similar (positive or negative) selection pressures, or if 
traits are positively correlated but submitted to antagonistic 
selection pressures. 

Recently, Agrawal & Stinchcombe [13] reviewed the impact of 
genetic correlations on the predicted rate of adaptation, gathering 
results from 45 studies on various plant and animal species, but 
found no general pattern: genetic correlations could either 
constrain or facilitate response to selection. However, two recent 
population-based studies that used the rate of adaptation metric 
defined by Agrawal & Stinchcombe ([13], see Methods) found that 
genetic covariances could decrease the rate of adaptation of life 
history traits by as much as 50% [14,15]. Hence, although it is 
difficult to generalize the impact of G on the response to selection, 
recent studies show that genetic constraints to evolution can be 
very strong. Presently, the scope for comparison of such metrics of 
multivariate constraints is very limited. As a result, more empirical 
work is needed across a range of species and traits in order to 
reach general conclusions about the influence of the G matrix on 
adaptation [16]. 

This study aims to generalize our knowledge of constraints or 
facilitation on the evolution of morphological traits, and to 
investigate how the interplay between selection and G matrices 
leads to such constraints. Our main objective is to estimate the 
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impact of genetic correlations on evolutionary trajectories using a 
comparative approach based on data from long-term field studies 
(> 12 years) of 10 populations of seven bird species. We chose four 
morphological traits (body mass, and length of tarsus, bill and 
wing, i.e., traits that represent shape and body size), because we 
were interested in assessing potential constraints for traits known to 
harbour considerable genetic variation e.g. [5,17,18]. First, we 
evaluate how G affects the predicted relative rate of adaptation, 
Ra [13] for these morphological characters. The rate of 
adaptation is the expected fitness gain due to an evolutionary 
response to current selection. The relative rate of adaptation (R A ) 
is obtained by comparing the rate of adaptation under models 
using observed G versus models with genetic correlations fixed to 
zero. Second, we seek to understand the origin of such patterns (1) 
by comparing the multivariate evolvability in the direction of the 
estimated directional selection (ep, which corresponds to the 
amount of predicted evolutionary response in the exact direction 
of current directional selection, P) versus the average evolvability 
in random directions of the phenotypic space, e [12,19]; and (2) by 
determining the orientation of the axis containing the highest 
percentage of additive genetic variance relative to the direction of 
selection (Fig. 1, [20]). This can be done by assessing the angle 
(9gmax) between g max and the directional selection p. These 
evaluations will allow the assessment of the extent to which a 
predicted micro-evolutionary response is facilitated or constrained 
by genetic correlations. 

Methods 

Ethical statement 

All data came from authorized monitoring of natural popula- 
tions and did not involve keeping birds in captivity. Such long- 
term studies require that birds are subject to minimal disturbance, 
and no manipulation was performed that would have caused 
animal suffering. Furthermore, all studies complied with national 
and international guidelines. All people collecting the data had 
banding permits. 

Species and focal traits 

Investigating evolutionary processes resulting from natural 
selection requires the use of data sets where phenotypes and 
relatedness are collected from populations in their natural 
environment. In this situation, estimating accurate G matrices 
necessitates long-term datasets with multigenerational pedigrees. 
We focused on four morphological traits that are most commonly 
measured in adult birds: wing length, tarsus length, body mass and 
bill length. Populations were hence selected based on the 
availability of a pedigree and the minimum number of morpho- 
logical traits needed. We gathered 10 data sets representing seven 
bird species from three continents (Table 1): red-billed gull 
(Chroicocephalus scopulinus, [21]), great reed warbler (Acrocephalus 
arundinaceus, [22]), barn swallow (Hirundo rustica, two populations, 
[23]), blue tit (Cyanistes caeruleus, three populations, [24]), collared 
flycatcher (Ficedula albicollis, [25]), Savannah sparrow (Passerculus 
sandwkhemis, [26]) and house sparrow (Passer domesticus, [27]). 

Wing length is a trait connected to flight performance and is 
especially important in migratory species [28] such as those 
included in this study (collared flycatcher, barn swallow, great reed 
warbler, Savannah sparrow). Tarsus length is a good approxima- 
tion for overall structural size in birds, because it is a skeletal 
measurement [29] . Body mass is also a general size measure, but 
more condition-dependent than tarsus length. Balbontin et al. [30] 
showed that body mass reflects condition, and that as such it 
provides a measure of changing condition among age classes and 
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generations. Bill length is associated with many characters, 
including foraging and song performance [31,32]. All of these 
traits have been shown to be heritable in several bird species [33- 
37]. 

Estimation of the additive genetic (co)variance matrix 

We estimated the G matrix in each population by using 
multivariate animal models [38,39]. Random effects included 
additive genetic effects (Va) and permanent environmental effects 
to account for repeated measurements of the same individual (V PE ) 
as well as a year effect (V year ). The analyses excluded measure- 
ments on offspring of the year. Age was included as a continuous 
variable (linear + quadratic) to account for aging effects on trait 
size. Tarsus length can change because of swelling or reduction of 
cartilage, wing feathers re-grow annually and are affected by aging 
and wear, beak length may become worn depending on diet, and 
body mass can be affected by age, e.g., because of decreased 
feeding performance. As we wanted to avoid losing power by 
removing individuals of unknown age (portion given in Table 1), 
we used mean substitution for individuals of unknown age: age was 
mean-centred and those individuals were assigned an age of zero. 
Because of power issues and technical complexity, males and 
females were not analysed separately so models contained sex as a 
fixed effect. When available and significant, we included a 
polynomial date effect (degree 2 or 3, according to significance) 
to control for mass and bill length variation during the breeding 
season. This affected the residual (co) variances, but not estimates 
of G. To avoid traits with larger means (Table 2) exerting a 
disproportionate effect on general patterns, we standardised traits 
prior to analysis. Because scaling to phenotypic variance (which 
can vary independently of additive genetic variance) can lead to 
problems of interpretation [40], we used standardization to the 
trait's overall mean [12,40]. 

A simple description of the multivariate animal model for one 
population is as follows: 

Y= n + Xb 4- Z a a + Z pe pe + Z yr yr + e ( 1 ) 

where Y is the vector of standardised phenotypic observations for 
all individuals, \i is a vector of mean phenotypes, b is the vector of 
fixed effects to be fitted (age, sex and date), and X is the design 
matrix relating phenotypic observations to the vector of fixed 
effects. Fixed effects were individually chosen for each population 
based on significance levels in a preliminary analysis (Table SI). 
For the random effects, a is the vector of additive genetic values, 
pe the vector of permanent environment effects, and yr the vector 
of year of measurement effect, with Z a , Z p( . and Z yr their respective 
design matrices. All random effects are assumed to be normally 
distributed, and elements of a are assumed to be drawn from 
a-A r (0,G®A) where G is the additive genetic variance-covariance 
matrix and A the relatedness matrix derived from the pedigree. 

All pedigrees were pruned using the R package "pedantics" [41] 
so they contained only informative individuals [41]. Details for 
each population are given in Table 1 and Fig. S 1 . 

Estimating selection 

To assess selection coefficients in each population, we used the 
classic approach by Lande & Arnold [42]. Directional selection 
gradients (P) were estimated by regressing relative fitness against 
morphological traits. Similarly, non-linear selection (y matrix) 
gradients were estimated using quadratic regressions, including 
cross products between traits, representing correlational selection 



PLOS ONE | www.plosone.org 



4 



Multivariate Constraints on Birds Evolution 



gradients. Quadratic coefficients from the regression were doubled 
so that they became analogous to selection coefficients [43] . 

Annual contribution to total individual fitness was estimated by 
yearly reproductive success (the number of fledged offspring). 
Morphological traits were first standardised by their means and 
then corrected for the same significant fixed effects as used in the 
animal models (i.e., effects of the fixed factors were subtracted 
from the actual measurement values), prior to selection analysis to 
obtain selection estimates consistent with the G matrices [12]. 

Each variable (fitness and morphological traits) was standard- 
ized within year, i.e. fitness was divided by annual population 
average success and we subtracted the mean annual phenotypic 
value from the overall mean standardized morphological variables. 

Estimating Constraints on and Facilitations of Responses 
to Current Selection 

We estimated how genetic correlations could affect both 
evolutionary trajectories and relative rate of adaptation in each 
population. First, we estimated the impact of genetic correlations 
on the predicted rate of adaptation in order to assess constraints on 
or facilitation of a response to the current selection acting in the 
populations, using the metric Ra defined by Agrawal & 
Stinchcombe [13]. This metric is the ratio between the predicted 
change in fitness given the predicted evolutionary change in mean 
phenotype per generation, in the presence of genetic correlations 
relative to what it would have been without these correlations. It is 
defined as: 



Ra 



Aw 0 (zy 



with 



A W(z) = Az r p + - Az T yAz . 



(2) 



(3) 



where A W{z) is the rate of adaptation (predicted change in fitness 
based on the predicted change of the mean phenotype of the 
population), Az the predicted change in average phenotype in the 
population calculated using the breeder's equation Az = GP, P the 
vector of directional selection gradients, and y the matrix of non- 
linear selection. In equation (2), AW g is the rate of adaptation 
taking into account genetic correlations, while A W 0 is the rate of 
adaptation when all the covariances between traits are set to 0. 
The ratio Ra is then compared to 1, with a ratio larger than 1 
implying higher rate of adaptation in the presence of genetic 
correlations (facilitation), while a ratio lower than 1 implies that 
genetic correlations slow down adaptation (constraint, [13]). To 
estimate the overall means across populations (equivalent to a 
meta-analytic mean) for relative rate of adaptation Ra, a ratio, we 
used the geometric mean: the overall mean was estimated for each 
iteration when estimating Ra for each population, so that a 
confidence interval could be built. 

Evolvabilities 

We also estimated multivariate evolvability and average 
evolvability [12]. Multivariate evolvability is the amount of 
predicted evolutionary response occurring in the exact direction 
of selection (ep, Fig. 1). It is estimated as 
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Table 6. Percentage of variance along g n 
the four morphological traits on g max . 



. and value of the first eigenvalue (xlOO) with 95% confidence intervals, and loading of 





Percentage of variance along g max 


First eigenvalue 


Wing 


Tarsus 


Mass 


Bill 


Red billed gull 


63.4 (58.28, 69.46) 


0.143 (0.118, 0.173) 


0.24 


0.415 


0.839 


0.257 


Great reed warbler 


53.66 (44.69, 64.6) 


0.123 (0.082, 0.166) 


0.117 


0.343 


0.669 


0.649 


Barn swallow - Badajoz 


62.97 (49.32, 72.15) 


0.159 (0.11, 0.23) 


0.178 


0.186 


0.898 


0.356 


Barn swallow - Kraghede 


63.38 (46.33, 74.2) 


0.222 (0.122, 0.335) 


0.251 


0.309 


0.916 


0.054 


Blue tit - Muro 


64.56 (48.99, 70.71) 


0.107 (0.069, 0.153) 


0.185 


0.388 


0.752 


0.499 


Blue tit - Pirio 


64.7 (55.52, 73.93) 


0.111 (0.073, 0.143) 


0.188 


0.474 


0.814 


0.279 


Blue tit - Rouviere 


59.38 (52.12, 65.7) 


0.127 (0.101, 0.16) 


0.174 


0.274 


0.864 


0.384 


Collared flycatcher 


59.5 (55.18, 64.45) 


0.131 (0.111, 0.151) 


0.154 


0.304 


0.919 


0.2 


Savannah sparrow 


62.24 (49.06, 68.75) 


0.163 (0.11, 0.217) 


0.145 


0.226 


0.8 


0.537 


House sparrow 


73.65 (60.72, 80.16) 


0.256 (0.171, 0.383) 


0.194 


0.494 


0.764 


0.366 



doi:1 0.1 371 /journal.pone.0090444.t006 



en = 



p r Gp 



(4) 



The angle between gmax and P cannot exceed 90° because 
gmax can be considered in its two opposite directions. Hence, if an 
angle larger than 90° was found, we took the complementary value 
180-6 mlax . 



Average evolvability over random selection gradients [12] 
represents the evolutionary potential associated with the G matrix 
if averaged across all possible directions in the phenotypic space. It 
is defined as 



e=E(X) 



(5) 



where Xs are the eigenvalues of G. Average evolvability thus does 
not depend on genetic correlations [12]. Note that our definitions 
of multivariate evolvabilities follow Hansen and Houle [12]. 
Evolvability can be denned as a univariate (variance scaled to the 
mean) or a multivariate estimate. Following [12,40], we use "Ia- 
evolvability" for univariate estimates of additive genetic variance 
scaled to the mean and "e" for multivariate estimates of 
evolvability. 

Angle between directional selection and g max 

gmax is the first eigenvector of G and the amount of additive 
genetic variance it contains is the eigenvalue of this vector. The 
sum of all eigenvalues of G represents the total additive genetic 
variance. Hence, the proportion of genetic variance along g max 
was estimated for each population by the ratio between the first 
eigenvalue of G and the sum of the four eigenvalues. This gives an 
assessment of the evenness of the distribution of the genetic 
variance in the different dimensions of G. 

The angle between gmax an <J the direction of selection (P, Fig. 1) 
estimates how close selection is from the axis that is the direction of 
least resistance. If selection and g max are aligned, the response to 
selection will be maximal while it will be constrained with 
increasing angles (with maximum constrain at 90°). The angle 
between g max and fi was calculated using: 



cosffL 



Smax P 
ilPllllgma, 



(6) 



Estimation method 

Both animal models and selection analyses were run using 
Bayesian methods with the MCMCglmm R Package [44]. The 
advantage of the Bayesian approach is that the use of posterior 
distributions facilitates the propagation of errors in estimates [45] . 
Although uncertainty around estimates of G matrices is usually 
large [8], attempts to integrate this uncertainty in the next steps 
(e.g., predicted response to selection) are extremely rare [14]. One 
of the goals of this analysis is to provide such estimates for each 
quantity described above. 

The posterior distribution was a sample of 1 000 values for each 
parameter. We used a total of 1,200,000 iterations for each 
analysis, with a burn-in phase of 200,000 and thinning of 1000. 
Priors were defined for variances and covariances. We assessed 
two priors for variances and covariances for each analysis: (1) a 
parameter expanded prior [46] and (2) a slightly informative prior 
(V = diag(n)*Vp/r, nu = n), where Vp is the phenotypic 
variance, n the number of traits and r the number of random 
factors. Our results were not sensitive to the choice of prior (Fig. 
S2). 

In the main text, we chose to present results from the model 
with slightly informative priors as it has a direct biological 
interpretation: the prior specification implies that (1) the variance 
is distributed evenly across the random terms and (2) traits are 
independent [44] . If information is coming from the prior, as it is 
built with null covariances, estimated genetic covariances would be 
biased downwards, if anything, and our estimates conservative. 

Results 

G-matrices 

iA-evolvabilities (lOOx additive genetic variance of traits scaled 
to the square of their mean), interpreted as the expected 
percentage of trait change per generation if it were submitted to 
selection as strong as on fitness itself [47], were on average 0.061 % 
(range 0.013%-0.178%, Tables 3, 4 and 5) across traits and 
populations. Heritability estimates were on average 0.30 (range 
0.05-0.60, Table S2). High I A -evolvabilities did not correspond to 
high heritabilities, and overall both estimates were unrelated 
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Figure 2. Relative rate of adaptation (R A ) in the 10 populations. Dots represent posterior mode estimates and lines the 95% confidence 
interval. The dotted black line at 0.75 represents the geometric mean of all populations while the dotted grey line at 1 shows the case of no effect of 
genetic correlations. Population number refers to the numbers given in Table 1. 
doi:10.1371/journal.pone.0090444.g002 



(R 2 = 0.04). The absence of congruence between evolutionary 
potential predicted from heritabilities and iA-evolvabilities is in 
line with a recent review [40]. 

Genetic covariances between all traits were positive in all 
populations (Tables 3, 4 and 5), and average genetic correlations 
were 0.35 (range: 0 to 0.76). In all populations, g max contained 
more than half of the total amount of additive genetic variance 
(geometric mean (95% CI): 61.3% (58, 64), Table 6) which 
suggests that G matrices were classically elliptical rather than 
spherical. The first eigenvalue, which represents maximal evolva- 
bility if selection and g max are aligned [12], was of the order of 0.1 
to 0.2 (values x 100, Table 6). All traits loaded positively on g max 
(Table 6), and body mass consistently had the highest loading on 
gmax< Because the first axis of a PCA can be interpreted as a size 
index, this suggests that the line of genetic least resistance (g m ax) is 
associated with body size. 

Natural selection on morphology 

Both the direction and strength of selection varied substantially 
across species but also across populations. In four of the 
populations (the three blue tit populations and the Kraghede 
population of barn swallows), directional selection was significant 
on bill length, tarsus or mass, but not on wing length (Tables 7, 8 
and 9). In collared flycatchers, great reed warblers and Savannah 
sparrows (Tables 7 and 9), we found significant directional 
selection on two traits, and in these three cases, selection was 
negative on mass and positive either on wing, tarsus or bill length, 
respectively. We also found no evidence of significant nonlinear 
selection. There was evidence for negative correlated selection on 
tarsus and mass in blue tits (Pirio) and barn swallows (Kraghede) 



and on tarsus and wing in barn swallows (Badajoz). Finally, there 
was significant positive correlated selection on bill length and wing 
in blue tits (Pirio), and wing length and mass in house sparrow. No 
significant selection was found in red-billed gulls. 

Constraints on predicted responses to current selection 

The predicted rate of adaptation was significantly lower in the 
presence than in the absence of genetic correlations (i.e., 95% of 
Ra values from the posterior distribution lower than 1) in four of 
the 10 populations (Table 10, Fig. 2): great reed warblers, blue tits 
in Pirio, collared flycatchers and Savannah sparrows. On average, 
Ra was 72%, which means that because of genetic correlations, 
the predicted fitness gain was on average 28% lower than it would 
be in the absence of these correlations. Despite large confidence 
intervals around the geometric mean across all populations, this 
average decrease was significant (geometric mean with 95% CI: 
0.72 (0.60, 0.85), Fig. 2), and no R A was larger than 1. 

Evolvability in the direction of P was on average 1.7 times lower 
than in random directions (mean ep xlOO (95% CI): 0.0369 
(0.0291; 0.0509), mean e x 100 (95% CI) : 0.0638 (0.0577; 0.0680), 
Table 10, Fig. 3), implying that current selection is acting in a 
direction of lower genetic variance than the average genetic 
variation in the phenotypic space. Confidence intervals within 
populations are much larger for ep than for e due to the 
uncertainty in the P estimates which adds to the uncertainty on G 
estimates. In accordance with these results on evolvability, the 
vectors gma* an d P were very close to orthogonal in most 
populations (Table 10), so that if all genetic variance was along 
ginax- no response to selection would be possible. However, other 
dimensions of phenotypic space include 40% of genetic variance, 
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so that multivariate evolvabilities (ep, the amount of response in 
the exact direction of P) were significantly different from zero 
(Table 10). These results emphasise that genetic variance 
remaining along dimensions other than g max also play a major 
role in determining evolvabilities. 

Discussion 

We report consistent evidence for multivariate constraints on 
morphological evolution across 10 avian populations studied in 
their natural habitat during extensive periods exceeding 1 2 years. 
Morphological traits generally display high heritabilities and 
harbour ample additive genetic variation [34,35,48-50]. There- 
fore they are often believed to be only weakly constrained in terms 
of evolutionary potential but see [40,51]. Here for linear 
measurements (mass excluded) we found iA-evolvabilities less than 
half (0.04% on average) of what was reported (0.09%) in the 
review by Hansen et al. [40]. The highest IA-evolvabilities were 
found for body mass, yet again for this trait, our estimates of Ia- 
evolvability were much lower (0.12 %) than the previously reported 
average of 0.94% [40]. In a univariate framework, for a trait with 
an iA-evolvability of 0.04%, this means that, if selection acting on 
this trait was as strong as on fitness itself, a change of 10% in the 
mean of the trait would be achieved in 240 generations [40] . 

Moreover, using a multivariate framework, we also found 
evidence of evolutionary constraints even when only four 
morphological traits were considered, emphasising that equating 
heritability with evolutionary potential can be misleading [40] . In 
fact, here we have shown that the predicted relative rate of 
adaptation (Ra) was on average 72%, which means that the 
predicted rate of adaptation was lowered by 28% (1- R A , range 
13-58%) due to the genetic correlations considered. 

Two scenarios may lead to a decreased rate of adaptation: 
negative genetic correlations with similar direction of selection 
pressures or positive genetic correlations in the presence of 
antagonistic selection. Negative genetic correlations have gained 
much interest in the study of evolutionary constraints [4] . This is 
mainly because selection is often positive on life history traits so 
that trade-offs should emerge as a consequence of negative genetic 
correlations for these traits but see [52]. However, genetic 
correlations between morphological traits generally seem to be 
positive ([53], this study, review in [54]). As the sign of selection on 
morphological traits is not always positive but depends on traits 
and populations ([49], this study), opposing selection patterns 
within the same organisms can be common and hence lead to 
constraints on responses to selection. Here, this scenario is 
illustrated by three populations of great reed warblers, collared 
flycatchers and Savannah sparrows, where the relative rate of 
adaptation was significantly lower than one. In these populations, 
antagonistic selection between mass and another trait (tarsus, wing 
and bill length, respectively), in the presence of positive genetic 
correlations explain this result. Such opposing selection patterns 
can arise because of selection for a specific function. For example, 
selection on wing length can be positive or negative, depending on 
whether long-distance flight or manoeuvrability are favoured e.g. 
[55]. Similarly, the sign of selection on beak size in Darwin's 
finches (Geospiza fortis) depends on the abundance of different seed 
types, which themselves depend on climatic events [56]. Further 
studies in each population would be needed to interpret selection 
patterns in terms of the function of traits, and to assess the 
ecological determinants behind these patterns. 

Such a reduction in the rate of adaptation reflects changes 
between the predicted responses to selection of traits whether or 
not genetic correlations are taken into account. In great reed 
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Figure 3. Comparison of evolvabilities in the direction of selection (e () , black symbols) and average evolvabilities in random 
directions of phenotypic space (e, grey symbols). Dotted lines represent the average value over the ten populations. Population number refers 
to the numbers given in Table 1. 
doi:10.1371/journal.pone.0090444.g003 



warblers, univariate models (i.e., not taking into account genetic 
correlations) predict significant responses in tarsus length and mass 
to selection, but no significant response in either trait is expected in 
the presence of genetic correlations (Table S3). In collared 
flycatchers univariate models predict a response to selection in 
both wing length and mass, but multivariate models predict a 
significant response only in mass. In Savannah sparrows, 
univariate models predict a response to selection in both mass 
and bill length, but only bill length is predicted to respond to 
selection in the presence of genetic correlations (Table S3). In 
contrast, no significant antagonistic selection was found in blue tits 
in Pirio, but multivariate models reveal nonetheless a disappear- 
ance of the response in mass when compared to univariate models. 
This is probably due to the fact that selection is significandy 
negative on mass while although non-significant, it is positive on 
the three other traits. 

We found a consistent pattern in that the orientation of g max 
was nearly orthogonal to directional selection in all populations. 
Although g max contained on average 60% of the additive genetic 
variance, the dimensions of G other than g max still contained ca. 
40% of additive genetic variance. This suggests that genetic 
correlations can decrease the rate of adaptation, but do not 
necessarily lead to an absolute constraint (i.e., here Ra^O). It is 
thus important to consider other dimensions along which additive 
genetic variance is distributed, and not only [19,57], as a 

reduction of the rate of adaptation of 28% is lower than what 
could have been expected based on the relative orientation of 
selection and gmax- 

In line with this argument, evolvability in the direction of 
selection (ep) was on average lower than evolvability in random 



directions of the phenotypic space (e), suggesting that selection 
may have reduced available genetic variance. This may be a very 
general pattern: depleted genetic variance in the direction of 
selection has also been found in sexually selected traits [20,58,59] 
and life history traits [60]. This result could suggest a depletion of 
additive genetic variance because of sustained directional selection 
on particular trait combinations [4]. However, there is still a 
debate about the stability of selection [61,62], so that a 
spatiotemporal analysis of selection patterns in each population 
would be required to assess whether sustained selection can be 
responsible for the observed pattern. Evolvabilities from this study 
(either ep or e) are very low compared to estimates from Simonsen 
and Stinchcombe [60] on life history traits of the ivyleaf morning 
glory (Ipomoea hederacea, ep~0.002) or foraging traits of three-spined 
sticklebacks (Gasterosteus aculeatus, ep~0.015 and e =0.007, [19]). 
However, results are similar to what was found by Bjorklund et al 
[63] in the same collared flycatcher population that we studied. 
There are still very few studies reporting estimates of multivariate 
evolvabilities, and it is not possible yet to interpret these differences 
either in terms of traits or taxa, yet we hope that our results will 
encourage further estimates in the near future. 

While our estimate of a decrease in predicted rate of adaptation 
(28%) is 2.5 times as large as the average decrease estimated by 
Agrawal & Stinchcombe [13] in their review (11%), they also 
found in 1 2 out of 45 studies that genetic correlations decreased 
the rate of adaptation by more than 30%. The decrease in 
predicted rate of adaptation from the present study is also lower 
than that found by Morrissey et al. [14] in their study of life history 
traits in a single island population of red deer (Cervus elaphus, 40%). 
In the Spanish population (Badajoz) of barn swallows, Teplitsky et 
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al. [15] found a decrease of 48% of the rate of adaptation for life 
history traits whereas in the present study we found a (non- 
significant) decrease of 20% for morphological traits for the same 
population. Two main factors may help explain such differences 
across studies. First, it is possible that morphological and life 
history traits differ in the amount of genetic constraints. Second, if 
selection is stable, constraints might actually be detected more 
readily than facilitation in natural populations. If genetic 
correlations facilitate the response to selection, populations should 
adapt and be subject to less intense selection. Hence, facilitation 
could be a transient state whereas constraints would represent a 
more stable state. Further analysis of data, such as those gathered 
in Agrawal & Stinchcombe [13], could provide valuable informa- 
tion as to when facilitation is more likely to occur. For example, 
does facilitation occur when organisms are subject to recent 
selection pressures, or when genetic architecture changes under 
new environmental conditions? 

The existence of multivariate constraints can have important 
implications for the potential of a micro-evolutionary response to 
rapid changes in the environment such as global climate change, 
because the pace of microevolution may be considerably reduced. 
The prevalence of such genetic constraints may begin to explain 
why so far little evidence of evolutionary adaptation to climate 
change has been reported [64]. The evolutionary significance of 
these constraints will also depend on the stability of the G matrix. 
The discussion regarding the extent to which and the conditions 
under which G is stable is still open, as some studies revealed 
either surprising constancy of G (review in [65,66]) or rapid 
changes [63]. 

Finally, our study showed significant multivariate constraints 
even though only four traits were included. This represents a very 
small fraction of all the traits integrated within an organism, and it 
is likely that constraints would become stronger if more traits were 
included [67]. As evidence is building that including more traits 
dramatically affects predicted responses to selection (e.g. [68], 
Table S3), and as our understanding and appreciation of 
evolutionary trajectories improves, it is becoming clear that 
multivariate studies should be the standard approach in evolu- 
tionary biology. Of course, including all traits is unachievable, but 
more comprehensive approaches based, for example, on modu- 
larity and identified suites of functionally related and highly 
correlated characters relatively independent of other suites of traits 
[69], promise to bring significant insights. 

Conclusions 

Our study assesses the generality of evolutionary constraints on 
morphology in birds that may arise from selection pressures such 
as those due to rapid environmental change. We found 
multivariate constraints on the predicted response to selection in 
morphological traits. Such traits are generally thought of as having 
a high evolutionary potential, which highlights the danger of 
equating heritability and evolutionary potential, as this can lead to 
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